x=c(1,3,2,3,3,3,4,6)
y=c(2,4,4,4,6,6,6,10)

X_bar = mean(x) # 辅助变量x的总体均值
y_lr_values = c() # 空数组，准备存放每个样本的回归估计

N=length(x) 
n=6 
CNn=choose(N,n)
AC=combn(seq(1,N),n) #所有样本的个体编号

for (k in 1:CNn){
xs=x[AC[,k]] #第k个样本：辅助变量
ys=y[AC[,k]] #第k个样本：主要变量

print(k)
print(xs)
print(ys)
x_bar = mean(xs) #样本均值
y_bar = mean(ys) #样本均值
b = cov(ys,xs)/cov(xs,xs) # 回归系数的最小二乘估计
y_lr = y_bar + b * (X_bar - x_bar) #计算回归估计
y_lr_values = c(y_lr_values, y_lr) #保存回归估计
}

print(y_lr_values)
print(mean(y_lr_values)) #打印所有回归估计的均值
print(mean(y)) #与上一个结果作比较
